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Abstract 

We study the structural properties of two-dimensional granular packings prepared by random 
deposition from a source line. We consider a class of random ballistic deposition models based on 
single-particle relaxation rules controlled by a critical angle, and we show that these local rules can 
be formulated as rolling friction in the framework of dynamic methods for the simulation of granular 
materials. We find that a packing prepared by random deposition models is generically unstable, 
and undergoes dynamic rearrangements. As a result, the dynamic method leads systematically to 
a higher solid fraction than the geometrical model for the same critical angle. We characterize the 
structure of the packings generated by both methods in terms of solid fraction, contact connectivity 
and anisotropy. Our analysis provides evidence for four packing regimes as a function of solid 
fraction, the mechanisms of packing growth being different in each regime. 
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I. INTRODUCTION 



Random ballistic deposition (RBD) is a well-known method for layer-by-layer construc- 
tion of random packings of hard particles such as granular beds and colloidal aggregates^ Q. 
This method is based on a simple and intuitive procedure. The particles (mainly monodis- 
perse spheres or disks) are allowed to fall sequentially along randomly positioned vertical 
lines over a horizontal substrate. Upon contact with the substrate or the first (already de- 
posited) particle, the particle either sticks or is further moved to a more favorable position 
according to a relaxation (or restructuring) rule. The RBD method can be efficiently imple- 
mented in a computer code for generating very large two- and three-dimensional packings. 
Elaborate large-scale simulations based on this approach have been used to investigate the 
geometrical properties of random packings (packing regimes, distribution functions, growth, 
etc)i0. 

It is obvious that the random deposition of particles can also be simulated by means of 
dynamic methods, such as molecular dynamics^, [f| and contact dynamics]!, [7j, §1, in the 
spirit of a real experiment where the grains are poured into a box. Such simulations require, 
however, substantially more computation time|9[]. This difficulty has been inhibiting enough 
to discourage for a long time systematic investigation of deposited beds following dynamic 
methods. But, the situation is far better today due to the fast increase of available computer 
power and memory during the last decade. There is now a considerable scope for dynamic 
simulations that can be exploited in order to study a number of highly interesting issues in 
the field of random packings. 

The objective of the present paper is to apply a contact dynamics algorithm to investigate 
randomly-deposited granular beds in a two-dimensional geometry. The geometrical texture 
(coordination number, solid fraction, ... ) of a granular bed depends on several physical 
parameters (particle properties, contact interactions, inertia of deposited particles) which 
can be tuned in a dynamic simulation in order to characterize the impact of each parameter 
on the texture. We propose here an original approach which allows to bring out some inter- 
esting features of granular beds in comparison to RBD models. We consider a generalized 
RBD model in which the relaxation of the falling particle upon contact with the substrate 
is controlled by the direction of the contact normal 9 (a single angle in 2D)|10[. Depending 
on whether 9 is below or above a critical angle 9 C , the particle either simply sticks or is 
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allowed to rotate until it reaches a local minimum position or forms a new contact below 6 C . 
Hereafter, we refer to this model as the CA (Critical Angle) model. The central feature of 



this model is that it allows to control the solid fraction p by varying the critical angle flit) 



The approach we propose consists in performing dynamic simulations of random ballistic 
deposition as closely as possible to the CA model. 

This requires that we transcribe the the above geometrical relaxation rule into a contact 
law which is reduced to the geometrical rule for the random deposition process. We first show 
that this requirement is met if the particles interact through a rolling friction law (similar to 
the Coulomb sliding friction) in which a contact torque is mobilized to resist relative rotation 
of two particles. We implement this law within a contact dynamics (CD) algorithm. Then, 
we perform two series of simulations both with the CA model and the CD method. In 
the first series, we use the granular beds prepared according to the CA procedure as initial 
configuration for a CD simulation. We show that, although local stability (sticking due to 
rolling friction or, alternatively, particles supported by two underlying contacts) is fulfilled 
for each particle added to the bed, the latter still undergoes collective rearrangements leading 
to a higher solid fraction. This implies that a granular packing prepared by geometrical rules 
is globally unstable. We study the extent of dynamic rearrangements and the structural 
properties of the CA and CD packings as a function of the critical angle. 

In a second series of simulations, we characterize the packings in terms of the average 
coordination number, structural anisotropy and contact connectivity as a function of the 
solid fraction. We show that the trends are globally similar for the CD method and the 
CA model (the same packing regimes are observed). We distinguish several packing regimes 
where different mechanisms (screening, chaining, branching, piling, jamming and ordering) 
are active and control the packing fraction. 



II. NUMERICAL APPROACH 
A. Critical-angle model 

Figure |l| shows the geometry of a contact formed by a falling particle i with a particle j of 
the substrate. The two particle centers define a line inclined at an angle 6^ to the vertical. 
For disks, % is also the direction of the contact normal n^, unit vector directed from the 



3 



center of particle j to the center of particle i. For brevity, we will refer to 6ij as the "contact 
direction" . 

The CA model is defined as follows ||T0f. If |^| is below a critical angle 8 C , defined in 
the range between (vertical direction) and 7r/2 (horizontal direction), particle i simply 
freezes by sticking to particle j. On the other hand, if \8ij\ exceeds 8 C , particle i is allowed 
to rotate around particle j until a second contact is formed with another particle k of the 
substrate. Then, there are three possible alternatives: 1) If \9ik\ < 9 C , particle i freezes as in 
the first case; 2) If the new position of particle i is a local minimum position, again particle 
i freezes; 3) If neither of the two latter conditions is fulfilled, particle i is again allowed to 
rotate around particle k until a new contact is formed with another particle of the bed and 
the three alternatives are examined again with this new contact. This procedure is iterated 
until particle i is stabilized either by sticking or by reaching a local minimum position. 

Two limits are of particular interest. When 6 C = 90°, all particles stick irreversibly 
to the substrate where ever they land. This limit corresponds to the random sequential 
adsoption model[14], [15|, [16], [HJ. When 6 C = 0°, all particles relax and the model is reduced 



to the steepest descent model [|l|, |12"| . In [13], the solid fraction was found to be p m in = 
0.3568 ± 0.0001 in the no-restructuring limit. The solid fraction for the steepest descent 
model is expected to come close to p max = 0.906 corresponding to a triangular packing. 
However, to achieve a structure with long-range ordering, the initial conditions are very 
important. In practice, the bottom line must initially be covered by an array of contiguous 
disks. Otherwise, simulations using the steepest descent algorithm have shown that the solid 
fraction will not exceed p = 0.82 which is the characteristic density of 2D monodisperse 
random close packing (RCP) where long-range order is broken by defects in the packing 
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In the CA model, the solid fraction p of the granular bed is a function of the critical angle 
6 C as shown in Fig. 0. In this figure, the solid fraction for each of the angles 9 C = 0, 1, 2, • • • ,90 
is an average over 30 independent runs. The bottom line was covered by a layer of 32 
contiguous particles and 1000 particles were deposited in each run. In order to avoid wall 
effects, periodic boundary conditions were implemented in the horizontal direction. Let 
us note that CA simulations are possible at much larger scales. Nevertheless, we restrict 
here the size and the number of CA simulations to those reasonably accessible to dynamic 
simulations since the results will be compared between these two methods below. Figure 
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shows that the solid fraction decreases monotonously from p m ax to p m in as 6 C is increased 
from to 90°. One can distinguish several regimes on this curve which will be discussed 
below in connection with dynamic simulations. 

The CA model is a geometrical model that meets the no-overlap condition between hard 
particles but involves a number of physical approximations about the stability of the packing 
and its growth. By nature, this model neglects inertia effects. The substrate is frozen and 
the relaxation step involves only the deposited particle. Moreover, the two stability criteria 
(local minimum position and sticking condition) for the deposited particle have a local nature. 
In other words, the model assumes that the whole packing remains in static equilibrium as 
long as all particles are sequentially stabilized by either of these conditions. 

In order to examine the validity of these assumptions, the approach followed in this paper 
is to perform dynamic simulations as closely as possible to the CA model and to compare 
the resulting packings. This implies that the particles should be released sequentially and 
they should hit the granular bed with negligibly small inertia. Moreover, upon contact with 
the substrate, the falling particle should dynamically stick or roll down, depending on the 
value of the contact angle with respect to the critical angle, until one of the two stability 
conditions is fulfilled. There is no difficulty in tuning the inertia in a dynamic simulation. 
But, we need to define a dynamic version of the relaxation rule. 

B. Rolling friction 

The dynamic content of stability due to a local minimum position is clear. The weight of a 
particle can obviously be balanced by the reaction forces exerted by two underlying particles; 
see Fig. 0(a). But the sticking condition requires both a contact force F^ and a "contact 
torque" Mij so as to counterbalance respectively the weight of the deposited particle 
and its moment with respect to the contact point. Figure |3|(b) illustrates this condition. 
The balance equations are 

Fij + niig = 0, (1) 
migriSmO + = 0, (2) 

where is the particle radius and g is the gravity. 

Let Nij and Ty be the components of the reaction force Fy along and perpendicular to 
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the contact normal n^. From equations § one gets 

Nij = rriig cos 9, (3) 
Tij = rriig sin0, (4) 
Mij = riN tj tan9. (5) 

The normal force Nj is positive (as it should) as long as —tt/2 < 6 < ir/2 (the angles 
are measured from the vertical). On the other hand, the relative sliding is inhibited if 
\Tij\/Nij < /j, s , where /x s is the coefficient of (sliding) friction (or equivalently, 9 < 9 S , where 
9 S = tan -1 ii s is the angle of friction). 

Now, if we require that particle i rolls only if 9 > 9 r , then from equations ^ we arrive at 
the following no-rolling condition: 



T"i Nij 



< tan6> r = fi r , (6) 



where /i r is a coefficient of rolling friction. Let us further assume that remains equal 
to its threshold value ±/i r rjA^j when rolling occurs. This condition is similar to the sliding 
condition = ±// s iVy. 

The rolling friction law, as defined here, and the more familiar Coulomb (sliding) friction 
law are shown in Fig. |] in the form of graphs |7|, |2(J. The rolling friction law relates the 
relative rotation velocity of the two particles to the contact torque My, whereas 

the sliding friction law relates the sliding velocity v*j to the tangential force T^. In fact, 
although for the sake of clarity we derived the condition ^| by considering the particular 
case of a deposited particle touching a particle of the bed, the application of rolling friction 
to a contact between two arbitrary particles in a packing is rather straightforward when 
formulated in the form of the graphs shown in Fig. |j. The torque transmitted through a 
contact to a particle in static equilibrium, for example, is the torque necessary to balance 
the sum of all force moments and other contact torques acting on the particle in the same 
way as the mobilized torque M in equations ^ counterbalances the moment rmg sin 9 of the 
particle weight. 

The prescription of rolling friction in a dynamic method follows the same steps as the 
sliding friction. The relation between the contact torque and the relative rotation velocity 
(Fig. |](a)) can not be represented as a motived function. Hence, in the framework of the 
molecular dynamics method, based on explicit integration of the equations of motion, this 



friction law has to be replaced by an approximate function^, |9], |2T]]. This "regularization" 
of the friction law is not necessary within the contact dynamics method which was employed 
for the present investigation |], [7], 

Using either of these dynamic methods, the sticking condition can be achieved if fi s is 
set to infinity (no sliding for no contact direction) and the angle of rolling friction 9 r is 
interpreted as the critical angle 9 C . This was implemented in our contact dynamics simu- 
lations. Alternatively, one may set /i r to infinity (no rolling for no contact direction) using 
the angle of sliding friction 9 S as the critical angle. It is also possible to use a combination 
of these two conditions. These conditions are not equivalent, but we will not discuss the 
differences in this paper. In all cases, the condition of sticking upon collision requires also a 
zero coefficient of restitution. 

To summarize, the following conditions allow to perform a dynamic simulation of random 
particle deposition in close analogy with the CA model: 

1. 9 r = 9 C 

2. No sliding (9 S = 90°) 

3. Weak inertia 

4. Zero coefficient of restitution 

The important difference is that, while in the CA model all degrees of freedom in the 
substrate are kinematically frozen, in our dynamic simulations only contact sliding is frozen 
by setting 9 S = 90°. All other degrees of freedom are active and the rolling friction governs 
all contacts: the contact between the deposited particle and the bed, as well as all contacts in 
the bed. Since the particles are not frozen in the granular bed, sequential particle deposition 
may thus lead to rearrangements in the granular bed. 



C. Simulation parameters 

The contact dynamics (CD) simulations involve a number of parameters which should be 
adjusted so as to minimize inertia effects without loosing numerical efficiency. The largest 
inertia are produced by the largest head-on velocity v max between colliding particles. Let 
At be the time step. The contact force due to inertia produced by a collision is mv max / At. 
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This force should be compared to the weight mg of one particle. Hence, the dimensionless 
parameter characterizing the ratio of inertia to weights is 

a = ^ (7) 
gAt 

The influence of a on the solid fraction and restructuring is an interesting subject in 
itself, but it will not be investigated in this paper. As emphasized previously, the focus here 
is put on those effects (equilibrium states, rearrangements) that arise from the geometrical 
configuration. Hence, we should use a low value of a. However, lower values of a mean 
slower simulations. Hopefully, the framework of the CD method allows for large time steps 
At up to the limitations related to the procedure of contact detection. On the other hand, 
the value of v max can be imposed for the falling particles, but further relaxation inside the 
packing may produce large impact velocities. In particular, at low solid fractions, where 
large voids exist in the bed, the free fall of a particle over distances compared to the system 
height H ~ 60r can give rise to impact forces far larger than the weight of a column of 
particles of the same height. 

In order to avoid such strong uncontrolled inertia, we implemented a "velocity barrier" 
trick that limits the particle velocities to v max = 0,3 ms -1 . With this choice, we can use a 
time step as large as At = 0.003 s. Then, setting g = 100 ms -2 , we get a — 1. This means 
that the largest impact force is just equal to the weight of a single particle. This choice is 
both reasonable and compatible with numerical efficiency 



III. DYNAMIC REARRANGEMENTS 
A. Stability of CA packings 

How stable are the granular beds prepared by means of the CA model? We have seen 
that the CD method, equipped with rolling friction together with suitable values of the 
parameters reducing inertia effects, meets the single-particle stability criteria of the CA 
model in the course of deposition, namely 1) the sticking condition as a function of the rolling 
friction angle 9 r (identified with the critical angle 9 C ), and 2) the local minimum position 
where the weight of a particle is balanced by the reaction forces at the two underlying 
contacts. Now, if we start a CD simulation using a packing constructed according to the 
CA model as the initial configuration, then one of the two following alternatives may occur. 
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If the single-particle stability criteria used in the course of deposition provide a sufficient 
condition for the global stability of the packing when the deposition is over, then the packing 
will remain in static equilibrium and the calculated forces will exactly balance all particles. 
Otherwise, the packing will be unstable and the CD simulations allow to calculate the 
particle rearrangements until a relaxed stable configuration is obtained. 

Our simulation data confirm rather the second alternative for nearly all values of the 
critical angle. One example is shown in figures jB|(a) and j5](b) for a packing of 500 particles 
with 8 C = 40°. The rearrangements occur in the whole bed, but they are much more hindered 
in the bulk than in the vicinity of the free surface. For this reason, the displacements appear 
mostly in the uppermost layers. The relaxed configuration has a larger packing friction. The 
solid fraction is still larger when the CD simulation is performed by random deposition of 
the same sequence of particles (as in the CA simulation) for the same value of 8 C (and the 
same boundary conditions), instead of using the CA configuration as the initial condition. 
The resulting packing is shown in Fig. |5|(c). In this latter case, the hindering effect related to 
the bulk density, which was active in the case (b), disappears since the CD rearrangements 
occur naturally in the course of deposition for each deposited particle. This means that the 
degree of instability of the CA packing shown in Fig. |5|(a) is more keenly reflected in the 
increase Ap of the solid fraction from (a) to (c) than from (a) to (b). 

Figure ^| displays the solid fraction p as a function of 8 C for packings prepared by CD 
sequential random deposition. The solid fraction for each value of 9 C is an average over 10 
independent CD runs. The CD simulations were performed for 40 different angels. The 
curve of p as a function of 8 C corresponds thus to 410 CD simulations of 1000 particles. 
Figure^ shows that, as expected, the solid fraction is every where larger for the CD method 
than for the CA model, except at 6 C = where a dynamic method requires an exceptionally 
high precision to reach a perfect triangular packing. Indeed, in this limit, tiny fluctuations in 
particle positions around a particle due to numerical overlaps are exponentially amplified in 
space as a result of long-range order |19[. Disregarding this pathological limit, the difference 
Ap is negligibly small for 9 C < 20° and 9 C > 80°. The largest variation Ap of the solid 
fraction, representing the largest dynamic rearrangements in the packing, occurs at 6 C ~ 50°, 
where p increases from 0.45 for the CA packing to 0.6 for the CD packing. 
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B. Influence on the packing structure 



The coordination number z (average number of contact particles around a particle), shown 
in Fig. as a function of 9 C , follows the same trends as the solid fraction. It is systematically 
larger in a CD packing than in the corresponding CA packing (for the same value of 9 C ) 
except in the 9 C = limit. A coordination number close to 2 reflects the predominance of 
particle "chains" in a highly porous packing. The coordination number increases from 2 to 
3 due to "branching" , and from 3 to 4 due to a growing interplay of chains. The increase of 



z beyond 4 requires long-range ordering! 22 1. This transition occurs only in the CA packing 
where the numerical precision is less stringent than in dynamic simulations. 

Due to dynamic rearrangements in CD random deposition, the CA and CD packings 
show also very different aspects as to the directional order of the contact network. The 
nonuniform distribution of contact directions can be characterized by means of the fabric 
tensor (j> defined from contact normals n k = (sin# fc , cosd k ) by|23|, |24[ 

1 Nc 

c k=i 

where N c is the total number of contacts, and n k a (resp. n%) is the a (resp. (3) component of 
the contact normal k. When the probability distribution function p(9) of contact directions 
is known, the fabric tensor can be calculated from the integral 

-tt/2 

n k a n k p(9) d9. (9) 

tt/2 



By construction, we have 0i + 02 — tr(<fi) = 1, where <f>\ and 02 are the eigenvalues. 
The mean contact direction in the packing is given by the major principal direction 9f of </>. 
The structural anisotropy of the packing is represented by a = 2(0i — fa)- The factor 2 is 
introduced in order to identify this value of a with that appearing naturally in a sinusoidal 
distribution p(6) = (1/tt){1 + a cos 2(9 - 6>/)}[gj|. 

Figure || shows the anisotropy of our granular beds as a function of 9 C . We see that 
the anisotropy of CD packings is systematically below that of CA packings. This is mainly 
because collective rearrangements tend to destroy columnar structures in a CD packing. In 
both cases, the anisotropy comes very close to zero for 9 C = 0. This effect is mainly related 
to the presence of a great number of particles with 5 and 6 contacts. In fact, using Eq. Rl it 
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can be shown that the anisotropy for the set of 6 contacts around a particle is zero, and for a 
set of 5 contacts around a particle can not exceed a threshold imposed by steric exclusions. 

The largest anisotropy in the CD packings is reached for 9 C = 90°, whereas the anisotropy 
of the CA packing passes through a maximum at 9 C ~ 50°. The anisotropy can be estimated 
analytically at 9 C = 90° where the packing growth is governed by sticking. Since the particles 
are released at random horizontal positions, the probability that a particle sticks at a contact 
angle 9 (with respect to the vertical) is p{9) = ^cos9. Note that the latter is a normalized 
probability density function over the range [— tt/2, 7r/2[. Using Eq. || with this expression 
for p(9), we find a{6 c = 90°) = 2/3 which is consistent with both CD and CA results at 
9 C = 90° shown in Fig. |J 

Since p{9) is an even function of 9 (p(0) = p(—6)), the major principal direction of the 
fabric tensor is vertical (6j = 0). However, this is only a consequence of symmetry and 
it does not imply that the distribution p{9) is peaked on 9 = 0. In fact, within each of 
the half-intervals [— tt/2,0] and [0,7r/2], the contacts have preferred directions. This can be 
seen in one example of p(6) for 9 C = shown in Fig. ^. We observe a local maximum at 
9 = 0, but there are local maxima also in each of the half-intervals. In order to extract the 
useful information about the direction of contacts, one can calculate the fabric tensor <fi by 
restricting the definition to one of the above two half-intervals. 

The major principal direction 9f for the interval [0, 7r/2] as a function of 9 C is displayed 
in Fig. [HI At 9 C = and 9 C = 90° both methods give the same direction, but every where 
else the contacts are more biased to the horizontal direction in CD simulations compared to 
CA simulations. This is an indication that the collective rearrangements reorganize contact 
directions. As for the anisotropy, the value of 6f can be calculated analytically in the 
limit 6 C = 90° over the interval [0, 7r/2] from the fabric tensor. The distribution function 
normalized over this interval is given by p(6) = cos 9 and the integral in Eq. |9| is calculated 
over the same interval. We find 9f(9 c = 90) ~ 32°, in agreement with the simulation result 
shown in Fig.[H| 

The existence of a local minimum in the evolution of 9/ for the CA model at 9 C ~ 30° 
or the changing of behavior in the CD curve at the same point can be understood as a 
consequence of competition between sticking and rolling. As 9 C is increased from zero, an 
increasing number of particles stick to the substrate at an angle in the interval [0, 9 C ]. These 
include both the ones sticking to the bed upon the first collision (whose number increases 
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as J Q c cos 8 d8 = sin 8 C ) and a number of the relaxed particles. On average, this subset of 
contacts tends to decrease 8f as long as 8 C is not too large. This explains the decrease of Of 
from 45° at 8 C = to ~ 28° at 8 C ~ 30°. But, the inclination of the contacts to the vertical 
increases at the same time following the increase of 9 C . This trend dominates clearly the 
evolution of Of beyond 8 C = 30°. 

IV. PACKING REGIMES 

The results presented in the last section show that, for a given value of the critical angle 
8 C , the solid fraction and the structure of the packing differ considerably from the CA model 
to the CD approach (excepted in the two limits of very loose and very dense packings). We 
attributed these differences to dynamic restructuring in CD packings as the particles are 
added to the substrate. However, in this section, we will show that the structure of a CA 
packing is quite similar to that of a CD packing if they are compared at the same solid 
fraction p (and thus, for different critical angles). This means that the structural properties 
of CA packings are quite realistic (close to CD packings) when they are considered as a 
function of the solid fraction rather than the critical angle. 

A. Fabric 

Figure [11] shows the coordination number z as a function of solid fraction p for CA and 
CD packings. In both cases, z increases with p. The two curves almost collapse for p < 0.6. 
For 0.6 < p < 0.8, the CD packings show only a slightly larger coordination number than 
the CA packings. For 0.8 < p, the CD packings show a slightly lower coordination number 
than the CA packings. The solid fraction p = 0.8 corresponds to z ~ 4 in both methods. 

The anisotropy a of the packings as a function of p is displayed in FigfT2]. The anisotropy 
decreases as a function of p for both methods except in the loosest CA packings where it 
increases a bit with p and passes through a peak before decreasing. The relatively low rate 
of decrease in the range p < 0.6 suggests that sticking is the dominant mechanism of growth 
in this regime, whereas rolling (or relaxation) is far more efficient in the subsequent range. 

Figure [13] shows the major principal direction 8/ of the fabric tensor restricted to the 
interval [0, 7r/2] (as defined in the last section) as a function of p. In CD packings, the 
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direction Of, representing the average direction of contact normals in the interval [0, 7r/2], 
increases quite slowly for p < 0.6 and much faster beyond 0.6. In CA packings, the anisotropy 
decreases for p < 0.6, passes through a minimum at p ~ 0.6 and increases in the subsequent 
range. However, let us remark that, as for the anisotropy, the difference in the value of Of 
between the two methods as a function of p is quite small as compared to the differences as 
a function of the critical angle (see Fig. [TDD • 



B. Connectivity 

The coordination number z is an average over all particles in a packing. But, the number 
of contact neighbors q varies in a packing from particle to particle. In a monodisperse pile, 
q can vary from 1 to 6. This "connectivity disorder" characterizes the disposition of the 
particles as "nodes" of the contact network. The connectivity of a packing is given by the 
fraction P q of particles having q contact neighbors. P\ corresponds to the "dead ends" 
of particle chains. The larger Pi, the stronger "screening" (the dead ends did not grow 
because they were screened by faster growing structures). P 2 and P 3 are related to chaining 
and branching, respectively. P4 corresponds to "piling", i.e. a natural situation where a 
particle is supported by two underlying particles and supports two others. P5 and Pq define 
"jammed" and ordered configurations. 

Figure |14] shows the connectivity numbers P q as a function of solid fraction for q varying 
from 1 to 6. The trends are globally similar in CA and CD packings and the differences for 
the two methods are quite small. All connectivity numbers vary monotonously with solid 
fraction except P3 which first increases to reach a maximum at p ~ 0.7 and decreases rapidly 
afterwards. In the range p < 0.6, P3 and P4 increase at the expense of Pi and P2 which 
decrease. P5 and Pq begin to increase significantly only at p ~ 0.7 and p ~ 0.8, respectively. 

The connectivity diagram P q is shown for three different solid fractions in Fig.|15|. The 
largest connectivity number is 2 for p < 0.6, 3 for 0.6 < p < 0.7, and 4 for 0.7 < p. 
Interestingly, the screening effect is more important in CD simulations (the CD curve for 
Pi stands above the corresponding CA curve). Chaining (P2), branching (P 3 ) and ordering 
(Pq) are slightly less important in CD simulations, while piling (P4) and jamming (P5) are 
enhanced. 

The above data show that the morphology of a CA packing is very close to that of a CD 
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packing at the same solid fraction. Both methods suggest four packing regimes characterized 
by the properties of the packing structure as a function of solid fraction: 

(a) p < 0.6: This regime corresponds to loose random packings characterized by chaining 
(i.e. P2 is the largest connectivity number), branching (P 3 increases as a function of p 
and becomes dominant at p = 0.6) and screening (P 1 is large). 

(b) 0.6 < p < 0.7: This is the regime of moderate random packings characterized by the 
largest value of P3 (branching) at the expense of chaining (P 2 ) which decreases rapidly 
as a function of p. 

(c) 0.7 < p < 0.8: This regime corresponds to dense close packings where piling is the 
main mechanism of growth and P4 is larger than other connectivity numbers. 



(d) 0.8 < p: This is the well-known dense ordered packing regime [18] characterized by 
z > 4. 

The transition densities (0.6, 0.7 and 0.8) appearing in this classification are approximate 
values. A more refined evaluation of these specific densities requires a considerably more 
computation time and a deeper insight into the mechanisms at play during the packing 
growth. 



V. CONCLUSION 



We investigated the structure of a class of randomly deposited granular packings whose 
density is controlled by a geometrical parameter, referred to as the critical angle. We used 
both a random ballistic deposition model with simple relaxation rules (the CA model), and 
a contact dynamics algorithm (the CD method) that incorporates those relaxation rules 
through a rolling friction law. The CD approach naturally leads to stable packings following 
dynamic rearrangements while in the CA model the packing is kinematically frozen after 
each single-particle relaxation. 

The following results were shown by means of extensive simulations: 

1) The packings prepared according to the CA model are generically unstable. When 
fed into the CD algorithm as initial configuration, the CA packings undergo dynamic 
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rearrangements. As a consequence, the solid fraction is larger in CD packings than in 
CA packings for the same critical angle (implemented as the angle of rolling friction 
in the framework of the contact dynamics method). 

2) The structural properties (anisotropy, connectivity) are quite comparable in CA and 
CD packings for the same solid fraction, even though significant differences were ob- 
served in packing anisotropies. 

3) Both methods reveal four packing regimes as a function of solid fraction, the prevailing 
mechanism of growth being different in each regime. 

An important outcome of this work is to show that the dynamic rearrangements are 
quite weak in the very loose and very dense limits where the structural properties are nearly 
the same. This means that, the random sequential adsorption model (irreversible sticking 
without relaxation), leading to very loose packings, and the steepest descent model (no 
sticking), leading to very dense packings, can be used with confidence in these two limits. 

The CD simulations reported in this work were meant to keep as close as possible to 
the CA model in order to perform comparable calculations with both methods. There is 
much more to be learned about the structure of particle packs generated by the CD method 
(or equivalently, molecular dynamics method). The influence of inertia and polydispersity 
on the observed packing regimes is currently under investigation. The shear resistance of 
deposited beds (e.g. in a biaxial compression) as a function of solid fraction is another 
important issue which we would like to address in near future. 
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FIG. 1: Geometry of a contact formed by a falling particle i with a particle j from the substrate. 
The contact angle 6ij is measured from the vertical. 
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(a) (b) 

FIG. 3: The two local stability conditions in the CA model: (a) local minimum position; (b) 
sticking. 
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FIG. 4: The graphs representing (a) sliding friction law and (b) rolling friction law; see text. 
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(c) 




FIG. 5: (a) An example of a CA packing; (b) The static packing obtained by the CD method 
starting with the CA packing in (a) as initial condition; (c) CD packing obtained by using the 
same sequence of falling particles as in (a). 
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FIG. 7: Coordination numbers z as a function of critical angle (in degrees). 
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FIG. 9: Polar diagram of the distribution p(6) of contact angles in a CD packing with 6 C = 0. 
The zero angle refers to the vertical direction. The angles are in degrees. 
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FIG. 10: The major principal directions 0/ of the fabric tensor as a function of critical angle. The 
angles are in degrees. 
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FIG. 11: Coordination numbers as a function of solid fraction p. 
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FIG. 12: Anisotropies a as a function of solid fraction. 
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FIG. 13: The major principal directions Of (in degrees) of the fabric tensor as a function of solid 
fraction. 
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FIG. 14: The connectivity numbers P±, ■ ■ ■ , Pq as a function of solid fraction. 
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FIG. 15: The connectivity diagram of CA and CD packings packings at three different values of 
solid fraction: (a) p = 0.5, (b) p = 0.7 and (c) p = 0.8. 
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